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ABSTRACT 

Stellar intensity interferometers will achieve stellar imaging with a tenth of a milli- 
arcsecond resolution in the optical band by taking advantage of the large light collect- 
ing area and broad range of inter-telescope distances offered by future gamma-ray Air 
Cherenkov Telescope (ACT) arrays. Up to now, studies characterizing the capabilities 
of intensity interferometers using ACTs have not accounted for realistic effects such 
as telescope mirror extension, detailed photodetector time response, excess noise, and 
night sky contamination. In this paper, we present the semi-classical quantum optics 
Monte-Carlo simulation we developed in order to investigate these experimental limi- 
tations. In order to validate the simulation algorithm, we compare our first results to 
models for sensitivity and signal degradation resulting from mirror extension, pulse 
shape, detector excess noise, and night sky contamination. 



1 INTRODUCTION AND GENERAL 
PRINCIPLES 

1.1 Intensity Interferometry 

Stellar Intensity Interferometry (SII) is an experimental 
method for measuring the angular diameters and acquir- 
ing high resolution images of stars. In a stellar intensity 
interferometer, the light from a star is received by two or 
more telescopes separated by a baseline that may range from 
tens of meters to kilometers, allowing for the resolution of 
surface features ranging from less than 0.1 milli-arcsecond 
(mas) (longest baseline) to 10 mas (shortest baseline) at 
visible wavelengths. Intensity interferometry relies on the 
correlation between th e light intensity fluctuations rec orded 
by different telescopes (|Hanburv Brown fc Tw iss 1957]). The 
fluctuations contain two components: the shot noise and 
the wave noise. The dominant component is the shot noise, 
which is the random fluctuation associated with the pho- 
ton statistics, and which is uncorrelated between telescopes. 
The smaller component is the wave noise, which can be 
interpreted as the beating between different Fourier com- 
ponents of the light reaching the different telescopes. The 
wave noise shows correlation between different telescopes 
provided there is some degree of mutual coherence in the 
light. 

The light intensity correlation between receivers 1 & 2 
is measured as the time-integrated product of the fluctu- 
ations 8i\ & 812 in the photodetector cu rrents. Althoug h 
higher order c o rrelations c a n be o f i nterest (iFontanal ( 
Gamd dl963) : ISato et alj i 19781 ); I Jain fc Ralston! 
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where (• • •) represents a ti me aver a ge. Ac c ording to the v an 
Cittert-Zernike theorem (|CittertJ l|l934 ): IZernikd l| 19381 ) ). 
the complex degree of coherence 7 is the normalized Fourier 
transform of the source radiance. One difficulty associated 
with image reconstruction using SII is that the measurable 
quantity is \^\ 2 , implying that the phase of the Fourier 
transform is lost. Methods of phase recovery and image 
reconst ruction have recently been developed and investi- 
gated (|Nunez et alj l|2012al ) and references therein) and 
Nunez showed that details of stellar surface features can 
be reconstructed w ith a relatively high degree of accuracy 
(|Nunez et al.ll2010l ). SII can provide stellar imaging to in- 
vestigate variou s topics of interest in cluding stellar rotation, 
limb darkening (I Nunez et al.ll2012b1 ), mass loss in Be-stars, 
surface temperature inhomogeneities (|Dravins et al.ll2010l ). 
and binary systems. 

It was shown by Robert Hanbury Brown and Richard 
Twiss that the signal S (photon 2 s -2 ) in an intensity inter- 
ferometry measurement is 



(iii • Si 2 ) = A 2 aV| 7 | 2 Ai/A/ 



(2) 



Ofir fc Ribakl (|2006al ). (|2006bl ). in this paper we will restrict 
ourselves to two-point correlations. Still, the Monte-Carlo 
simulation approach described here could also be directly 
used for higher order correlation studies. 

In the case of a thermal light source and an ideal in- 
tensity interferometer, the two-point correlation is equal to 



where A is the effective light collecting area of the tele- 
scopes, a is the quantum efficiency of the photodetectors, 
r\ (photon ■ m -2 • s _1 • Hz -1 ) is the spectral density of the 
light, A/ is the signal bandwidth of the photodetectors and 
electronics, Av is the o ptical bandwidth of th e light, and T 
is the integration time l|Hanburv Browr]|l974l ). The noise N 
(photon 2 s~ 2 ) is 



N = V2Aa V Av(Af/T) 1/2 . 



(3) 
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Therefore the signal-to-noise ratio SNR for an intensity 
correlation measurement with an ideal system is 



S 



SNR=-=Aar 1 \~ i \ 2 y/AfT/2 



(4) 



It is worth noting that the SNR is independent of the optical 
bandwidth Av used for the observations. This expression for 
the SNR only accounts for the photon statistics, previously 
referred to as the shot noise. The factor of v2 accounts for 
the fact that starlight is non-polarized. If we consider fully 
polarized light then the SNR is increased by a factor of v2- 

There is increasing interest in SII because of the rela- 
tive ease to achieve long baselines as well as the advantages 
it may offer in terms of cost if implemented in conjunction 
with existing or futur e very h igh energy (VHE) ga mma-ray 
ACT arrays (|Weekesl l|2003h : iDravins et all (|2012l )). ACTs 
operate by taking advantage of the Cherenkov light flashes 
produced by atmospheric showers. Because the atmospheric 
Cherenkov light is very faint, observations of VHE gamma- 
rays require large light collectors (A ~ 100 m 2 ) and are re- 
stricted to nights with little or no moonlight. During moonlit 
nights, when ACT arrays are less effective for gamma-ray ob- 
servations, they could be used for SII observations through 
narrow optical bandwidths. The implication is that SII mea- 
surements could be performed with ACT arrays with mini- 
mal interference with the VHE observation programs while 
increasing the scientific output of these instruments. 

Arrays such as the Cherenkov Telescope A rray (CTA), 
which could consist of up to ~100 telescopes (|Actis et al.l 
2011) could provide thousands of different baselines simulta- 
neously and could offer detailed imaging capabilities. In pre- 
vious studies characterizing SII sensitivity and imaging ca- 
pabilities, it has been assumed that the detectors are point- 
like in size so the van Cittert-Zernike theorem applies di- 
rectly and the degree of correlation between the intensity 
fluctuations is strictly proportional to the squared magni- 
tude of the Fourier transform of the source radiance. How- 
ever, when the telescope aperture becomes comparable to 
the baseline separating the telescopes, then the light may 
no longer be regarded as fully c oherent across individual 
apertures (tHanburv Browr]|l974l 'l and the correlation data 
then departs from a pure Fourier transform. 

Additionally, other instrumentation-related systematic 
effects have been neglected in previous studies. For exam- 
ple, electronic artifacts such as single photon response pulse 
profile and excess noise affect the signal or contribute to 
the degradation of the SNR. Also, the effect of the night 
sky background light integrated in the point spread func- 
tion (PSF) of the light collector and the profile of the op- 
tical band pass have only bee n qualitatively or very ap - 
proximately taken into account l|LeBohec and Ho lder 2006). 
Furthermore, high speed digitization electronic systems are 
considered to be used for SII applications as they offer the 
flexibility of offline signal correlation analysis. The develop- 
ment of data analysis algorithms for such systems requires 
realistic simulated data. Finally, other effects such as mirror 
non-isochronism and inaccuracies in time delay lines and 
star tracking can be investigated. All of these aspects are 
important for the performance characterization, data anal- 
ysis preparation, and the design and deployment of an SII 
observatory. 

In this paper we present a Monte-Carlo simulation of 



a semi-classical quantum description of light and a simple 
instrumentation model we developed in order to investi- 
gate the above mentioned instrumentation-related effects. 
Section [2] describes the simulation approach to achieve the 
proper photon statistics at the different telescopes and it also 
includes models of instrumentation effects such as the single 
photon response pulse and excess noise. Section[3]presents a 
few simulation applications characterizing instrumentation 
related effects, concentrating on the finite telescope diam- 
eter, the single photon pulse shape, the excess noise, and 
the night sky contamination. Since there is no actual data 
available to test our simulations against, in each case we 
compare the results to simple models corresponding to ideal 
cases to obtain further validation of our approach. Finally, 
the findings are summarized in Section [4] 



2 PHOTON LEVEL MONTE-CARLO 
SIMULATION OF SII 

2.1 Thermal light source and telescope signals 

We model a stellar light source with a total photon flux 
$* (photon ■ m -2 • s _1 ) within a given optical bandwidth 
Av as a collection of M discrete sources, each contributing 
an equal flux $a/ = $*/M. Each point source is defined 
by a wave amplitude Aj , an angular frequency uij , a phase 
4>i, and an angular position 0j — (8 x ,8 y ). At a given time, 
the phase <j>j of the light is taken randomly and uniformly 
from [0, 27r] and t he amplitude A 1 is tak en randomly from 
a Gauss deviate (|Mandel fc Wolj|l995l ) of mean ^7 = 
and variance _4 2 = $m/c, where c is the speed of light. 
In order to make the simulation closer to the continuous 
distribution of a realistic stellar source, at each time-step 
we also randomly set the angular frequency luj and angular- 
position 8j — (6 X , 9 y ) of each point source from distributions 
corresponding to the spectral density spectrum, radiance, 
and size of the simulated stellar source. 

Each telescope mirror is modeled as a set of small light 
collecting elements of area dAt centered on position Xu = 
(xk,Uk)- Note that "small" in this context means that the 
mutual degree of coherence between any two points within 
any single area element is maximal. 

The average number of photons emitted by the star and 
incident on one area element k of the telescope during a 
time-step St can be written: 



d/it, = dAk ■ St ■ c 



l )+<t>j 



(5) 



Throughout this paper, the notation denotes a statistical 
average while (■ • ■) denotes a time average. 

The average number of star photons collected by tele- 
scope i during a given time-step can then be written: 



(6) 



Note that (fk) is the time average of fii and (p,i) = n i , i where 
n- ki — Ai St is the average number of photons emitted by 
the star and collected by the telescope in one time-step. 
The summation of the waves prior to taking the magnitude 
is responsible for the non-Poisson distributed light intensity 
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fluctuations, previously referred to as wave noise, displaying 
correlation between the different telescopes. 

The actual number of photons reaching telescope i in 
a time St consists of light from the star as well as stray 
background light and can be written as 



IlT, = P(/M + UBGi) = P(jM + P n i , i ) 



(7) 



where P(m) represents a Poisson deviate of mean and vari- 
ance m and the night sky contamination is defined in terms 
of the average source radiance (tibg^ = Pn* i where /3 is a 
dimensionless factor which sets the amount of stray light). 

In principle, with the description provided thus far, 
we can simulate any intensity interferometry signal, tak- 
ing the time-step to be smaller than the coherence time 
t c (~ 1CP 5 ns) of the light. However, in practice, it is de- 
sirable to be able to use a time-step St ^> t c comparable to 
the electronic time resolution 5t e (~ 1 ns) while maintaining 
a manageable computation time. 

To do this, the correlated, non-Poisson distributed pho- 
tons are artificially diluted with a stream of purely Pois- 
son distributed photons. Contaminating the starlight with 
purely Poisson light degrades the SNR equivalently to de- 
grading the signal bandwidth, permitting simulations with 
5t^> t c . Without affecting the mean number of photons inci- 
dent on the telescope, the instantaneous number of photons 
may be written: 



71t ; = P(n* 4 (l — k) + jUj/t + /3n* 



(8) 



where the parameter k sets the degree of dilution of the 
correlated photons without affecting the total photon rate. 
k — corresponds to a case in which the correlation is zero, 
while k = 1 corresponds to a case of maximal correlation. 

We want to obtain the correlation between the signal 
fluctuations about the mean. They can be written = 
riT i — UTf, where we use fiT i = + /?) so as to ensure 

(»T < )=0. 

The correlation is computed as 



g2 = & 



S T 2 ) 



(9) 



After correction of the mean of g 2 for the parameters k & /3, 
(see Appendix f5J, we see that 



(10) 



in the limit where the diameter of the telescopes is small 
compared to the distance required to begin to resolve surface 
features of the star. Note that k = y/TcJSt is equivalent to 
the y/AfT term in Eq.g]with A/ = l/8t e . 

Photons are recorded using electronic systems with a 
specific time response. It should be noted that in SII, it 
is the correlation between the high frequency fluctuations 
in light intensity recorded by different telescopes that is 
measured. Such high frequency fluctuations are directly ac- 
cessible through an AC coupling of the photo-detector. An 
AC coupled trace can be modeled using any desired single- 
photon response pulse whose integral is zero. The accumula- 
tion of individual photons, each modeled by an AC coupled 
pulse, directly provides the telescope signal st ( . With a de- 
tailed pulse model, the sensitivity of a given system may 
be evaluated more precisely and the effects of experimental 
timing inaccuracies may be investigated. In this paper, for 
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Figure 1. The simulated data for the observation of a uniform 
disk star, 2 mas in diameter is shown to reproduce an Airy disk 
profile. See text for details. The data sets shown here are: a star 
consisting of 10 point sources which are not randomized at each 
time-step (these points are indicated by □); a star consisting of 
10 point sources which are randomized at each time-step (indi- 
cated by *); and a star consisting of 100 point sources which are 
randomized at each time-step (indicated by +). Deviations at 
large baselines are greatest when the star is simulated with fewer 
points and when the points of the star are not randomized at each 
time-step. 



the application of our simulations, we restrict ourselves to 
studying the simplest case of an AC coupled square pulse 
so that our results can be compared to simple models to 
validate our approach. 



2.2 Proof of principle of the light model 

To test that the simulations yield the correct correlation, 
we simulate a uniform disk star and calculate |7| 2 as a 
function of the baseline for point-like telescopes, (using a 
square pulse of width 8t e ), and neglecting detector excess 
noise and stray light contamination. Figure [1] shows the 
simulated data for a uniform disk star, 2 mas in diameter 
observed with a point-like telescope with a 1 m 2 flux col- 
lection area through a 10 nm optical bandwidth centered 
on a wavelength of 400 nm. The star, which has a flux 
= 10 9 photons • m~ 2 • s _1 through the selected optical 
bandwidth {my ~ —0.5), was simulated for a duration of 
1 ms. The unrealistically large signal-to-noise ratio is ob- 
tained because of the use of a signal bandwidth that is only 
10 times smaller than the optical bandwidth. We see that the 
signal reproduces an Airy disk profile to high precision. How- 
ever, we found that the fewer the number of point sources 
making up the star, the greater the power at high frequen- 
cies, which causes |7| 2 to deviate from the Airy disk profile. 
Additionally, we found that randomizing the location of each 
point source within the spatial extension of the star as well 
as randomizing the frequency emitted by each point source 
at each time-step is more representative of a realistic star 
and further reduces deviations from the Airy disk. 

To further test our approach, we compared the expected 
SNR calculated from Eq. [4] to the SNR of simulated cor- 
relations. In these experiments, the excess noise and stray 
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Figure 2. Statistics of simulation results are shown to follow the 
theoretical signal-to-noise prediction by Hanbury Brown. The + 
shows data acquired from varying the baseline, the X shows data 
acquired from varying the pulse width St e , and the * shows data 
acquired from varying the flux The ratio of the simulated 
to expected SNR for various parameters fall around a 1:1 ratio 
indicated by the straight line. 



light contamination are set to zero. Each numerical exper- 
iment was run 50 times in order to calculate the standard 
deviation of each set of experimental parameters for a non- 
resolved star (i.e. 7I 2 — 1). In Figure [2] we compare the 
SNR obtained from simulations to the SNR calculated from 
Eq. [4] For the signal bandwidth, a square pulse was used so 
that the width of the pulse is unambiguous. The AC cou- 
pling was taken into account by subtracting the average sig- 
nal rather than including a negative tail for the pulse, so the 
correlation is obtained between signals st ; as described in 
Section \2. II The width St e of the pulse was varied from Ins 
to 30 ns and it is verified that the standard deviation evolved 
as a square root of the pulse width. Similarly, the flux of the 
light source was varied from $* = 10 8 photons • m -2 • s _1 
to = 10 9 photons • m -2 ■ s _1 . We see that the simulated 
SNR precisely follows the prescription of Eq. [4] The effect 
on the SNR of varying the flux is equivalent to varying the 
light collecting area of the detectors. The effects of increas- 
ing the aperture will be further discussed in Section 13.11 
These results confirm that our Monte-Carlo approach pro- 
vides consistent results in simple cases, and so, more subtle 
instrumentation artifacts can be investigated. 



3 SIMULATIONS TOWARD A NON-IDEAL 
INTERFEROMETER 

3.1 Mirror extension 

When the effect of the telescope mirror extension is impor- 
tant, corrections must be applied to the data before the van 
Cittert-Zernike Theorem can be used. Simulations can be 
exploited to establish this correction. The mirror extension 
effects are significant and corrections must be applied in 
most cases since telescopes that may be constructed in fu- 
ture arrays such as CTA could be up to ~ 30 m in diame- 



ter l|Actis et al.ll201lh , which is comparable to intertelescope 
distances of current and future arrays. 

We have investigated the effect of mirror extension on 
a 2 mas diameter uniform disk star for various telescope di- 
ameters which are comparable to the dish sizes considered 
by CTA. We see that the shape of the ideal curve, which is 
shown by the Airy disk profile and by the simulated data 
for the case of a point-like detector (i.e. too small to start 
resolving the star), is smoothed out as the size of individual 
detectors increases and begins to resolve the star. 

In principle, the effect of large detector sizes on the de- 
gree of correlation is equivalent to taking a successive con- 
volution of |7| 2 with the shape of the light collecting area of 
each individual telescope (see Appendix[B} ( Hanbury Brown! 
Il974h . which moves \y\ 2 away from being the squared mag- 
nitude of the Fourier transform of the source. Because the 
simulated star is a uniform disk, \j\ 2 follows an Airy disk 
profile independently from the orientation of the direction of 
the baseline with respect to the star. Therefore, the Airy disk 
in two dimensions may be obtained by assuming axial sym- 
metry in the correlation data. In order to test the simulation 
algorithm against the successive convolution model, we sim- 
ulated a pair of identical telescopes. The double convolutions 
of the simulated two-dimensional data for two telescopes 
with uniform disk-shaped light collecting areas of equal di- 
ameters, 10 m, 20 m and 30 m respectively are shown in Fig- 
ure Within the standard error, each data set agrees well 
with the prediction. 

The effect of mirror extension on image reconstruction 
capabilities requires a detailed study which cannot be car- 
ried out in this paper. It may be possible to develop correc- 
tion algorithms to apply to the data before analysis to par- 
tially alleviate the effect; however, some information is lost 
in the successive convolutions. Alternatively, when a param- 
eterized image model is available, detailed simulations can 
be compared to the data for parameter optimization. 

3.2 Excess noise 

In the high-frequency regime of interest here, most of the 
fluctuations of the telescope signal correspond to the Pois- 
son statistics of the collected photons. To this, we must add 
the fact that the detector response may fluctuate from pho- 
ton to photon. This is known as the excess noise. For exam- 
ple, in the case of a photo-multiplier, the excess noise results 
primarily from the fluctuations in the number of electrons 
ejected at the first dynode. A typical excess noise level for a 
photo-multiplier is around 30% of the average single-photon 
response. The excess noise is uncorrelated between differ- 
ent telescopes but it may have the effect of reducing the 
SNR achieved in the detection of a correlation between tele- 
scopes. Here, we model the excess noise by multiplying each 
single photon pulse by a Gaussian deviate whose mean is 1 
and whose standard deviation is a. However the Gaussian is 
truncated at zero to avoid multiplication of the pulse by a 
negative factor. As a consequence, the mean is greater than 
the most probable amplitude. In order for the mean single 
photon response amplitude to remain the same, all signals 
are divided by the mean of the excess noise distribution. 
The SNR expression in Eq. U accounts only for the fluctu- 
ations associated with photon statistics, so the simulations 
are used to gain an understanding of the noise introduced by 
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Figure 3. An Airy disk profile is shown by the solid line. The 
double convolutions of the Airy disk profile with uniform disks 
of diameters 10 m, 20 m, and 30 m are shown by the dashed 
lines. The double convoluted profiles are plotted with simulated 
correlation data for telescopes with mirror diameters of 10 m (x), 
20 m (*), and 30 m (□) for comparison. Aside from the telescope 
extension, the physical parameters used in the simulations arc the 
same as used in Figure [l] 




Excess noise 

Figure 4. Dependence of SNR on excess noise. Simulations were 
made using a pair of 1 m 2 telescopes and an observation time 
of 1ms. For each point, the SNR was obtained from running the 
simulation 500 times. The + indicates data for fluxes of = 
10 6 photons • m~ 2 ■ s — 1 (approximately corresponding to a visual 
magnitude my K. 7 through a lOnm optical bandwidth), * indi- 
cates data for fluxes of <!>* = 10 7 photons ■ m~ 2 ■ s _1 (my 4.5), 
and X indicates data for fluxes of = 10 s photons ■ m — 2 ■ s _1 
(my RJ 2) and the corresponding solid and dashed lines indicate 
the respective models. The model of the SNR fits the simulated 
data well for lower values of excess noise. 



the electronics. Figure [4] shows the SNR obtained from the 
simulations and from a model developed in Appendix [Cl as 
a function of the excess noise. The model does not account 
for the truncation of the distribution of the single photon 
response pulse at zero. Therefore, as expected, simulated 
data deviates increasingly from the model as the relative ex- 
cess noise is increased. The simulated SNR does not degrade 
as fast as suggested by the model since the distribution of 
pulse amplitudes is narrower due to the truncation effect. 
Since the difference in SNR between simulations with = 
10 6 photons • m~ 2 ■ s" 1 and = 10 7 photons ■ m~ 2 ■ s _1 at 
zero excess noise is smaller than the statistical error, we in- 
creased the number of simulation experiments by a factor 
of ten while maintaining all the parameters as in Figure [3] 
in order to establish a clear difference between the sensitiv- 
ities for these fluxes. We found that from the simulation, 
the normalized SNR for a flux = 10 6 photons • m -2 • s _1 
is S/N = 0.706 ± 0.007 compared to a model prediction of 
S/N = 0.706 and for a flux = 10 7 photons ■ m" 2 ■ s" 1 , 
the SNR is S/N = 0.685 ± 0.007 compared to a prediction 
of S/N = 0.693. For non-zero excess noise, the SNR departs 
from the evaluation using Eq. [4] as discussed above. This 
type of study allows for a more accurate and realistic esti- 
mate of the sensitivity of an SII experiment. A more realistic 

model of the excess noise statistics, such as a Polya distri- 

a II 1 

bution for photo-multiplier tubes (|Prescottlll966n can easily 
be incorporated into the simulation algorithm. Similarly, the 
noise characteristics of other types of detectors, such as sili- 
con photo-multipliers (SiPM), Geiger-mode avalanche pho- 
todiodes, and micro- channel plate (MCP) photo-multipliers 
i|Wagner et al.ll2009Ti could be modeled as well. 



3.3 Night sky contamination 

The simulations can also be used to investigate the effects 
of stray light contamination. When very faint stars are to 
be observed, the observation time must be long enough to 
obtain a sufficient SNR for the measurement. In these cases, 
contamination from the night sky background, especially 
during bright moonlit nights, cannot be ignored. When the 
light received is dominated by the night sky background 
contamination, increasing the observation time is no longer 
beneficial. An arbitrary amount of stray light can be in- 
cluded in the simulation as a stream of additional photons 
with purely Poisson statistics. The expected dependence of 
the SNR on the amount of stray light is derived in Ap- 
pendix [X] To test this, we simulated a star of brightness 
= 10 8 photons • m~ 2 ■ s _1 and integrated over T = 0.1 /js. 
The dependence of the SNR as a function of the stray light, 
in terms of the brightness of the source, is shown in Figure 
[S] We observe a good match between the simulated data and 
our model. 

Stray light contamination ultimately limits the faint- 
ness of stars that can be observed. Using the brightness 
of the night sky for various phases and distances from the 
moon from iKrisciunas Sz Schaeferl |l99 ll ). Figure [6] shows 
the magnitude of the stars whose brightness equals the 
night sky luminosity integrated within the PSF. This sit- 
uation corresponds to j3 = 1 which results in a degra- 
dation of the SNR by a factor of 2 compared to ft — 
(i.e. no background light). We used an optical PSF full 
width at half maximum (FWHM) of 0.06°, correspo nding 
to the VERITAS telescope array (|Holder et all 120061 '). An- 
other limitation for the observability of a star is the prac- 
ticality of the required integration time. For this reason, 
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Figure 5. Dependence of SNR on the amount of background light 
f3 in terms of the source flux = 10 s photons ■ m~ 2 ■ s —1 from 
an unresolved source observed for a duration of T = 0.1 /is with 
1 m 2 light collecting area telescopes. For each point the SNR was 
obtained from running the simulation 500 times. 



in Figure [6] we include the integration time required to 
achieve a SNR = 5 for a star with spectral density rjv = 
(5 ■ 1(T 5 photons ■ m~ 2 • s" 1 ■ Hz" 1 ) • 2.5~ v at |7| 2 = 0.5 us- 
ing a telescope of area A — 100 m 2 , a photodetector quan- 
tum efficiency of a = 0.25, and a signal bandwidth of 
A/ = 100 MHz. This does not take into account the fact 
that in a telescope array, the redundancy of the baselines 
can be used to improve the sensitivity. 

In addition to limitations imposed by the moon, SII 
measurements may also be affected by other stars in the 
field of view of the telescope, especially when observing very 
faint stars. If there is a star which is as bright or brighter 
than the star of interest in the PSF of the telescope then it 
is no longer possible to measure the star of interest. As the 
brightness of the star decreases then the average number of 
stars within the PSF of the telescope increases, as shown 
in Figure [7] However, for stars brighter than magnitude 12, 
given the low density of bright stars, other stars within the 
telescope PSF are unlikely to interfere with measurements. 
Stars of magnitude 12 or fainter will likely remain out of 
reach of SII because of the required integration time. 



4 CONCLUSION 

Previous studies have not accounted for many realistic 
effects associated with SII measurements. Here, we have 
demonstrated that our simulations correctly reproduce our 
simple models predicting how \ j\ 2 and the measurement sen- 
sitivities evolve with various parameters, and we have briefly 
discussed instrumental properties such as telescope mirror 
extension, signal bandwidth limitation and electronic pulse 
shape, excess noise, and stray light contamination. We have 
verified that the signal measured with spatially extended 
detectors is the successive convolutions of the normalized 
Fourier Transform of the source with the shape of each de- 
tector area. We also verified that the SNR degrades as the 
excess noise increases according to a simple model devel- 
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Figure 6. Faintest stars that can be observed using SII with 
telescopes with a PSF of 0.06° FWHM under various moonlight 
conditions. The horizontal axis shows the angular distance of the 
star from the moon, the vertical axis shows the faintest visual 
magnitude star that can be observed, and the various curves cor- 
respond to different phases of the moon, where + indicates a 
moon phase of 30 deg, X indicates a moon phase of 60 deg, * in- 
dicates a moon phase of 90 deg, and □ indicates a moon phase 
of 120 deg. The right vertical axis shows the integration time re- 
quired to achieve a SNR = 5 for a star with spectral density 
nv = (5 ■ 10" 5 photons • m~ 2 ■ s" 1 • Hz" 1 ) • 2.5~ v at \-y\ 2 = 0.5 
using a telescope of area A = 100 m 2 , a photodetector quantum 
efficiency of a = 0.25, and a signal bandwidth of A/ = 100 MHz. 
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Figure 7. The vertical axis shows the average number of stars in 
the field of view of a telescope with a PSF FWHM of 0.06° as a 
function of the apparent visual magnitude V. The various curves 
correspond to the angular distances from the galactic plane, where 
+ corresponds to an angular distance of Odeg, X corresponds 
to an angular distance of 30 deg, * corresponds to an angular 
distance of 60 deg, and □ corresponds to an angular distance of 
90 deg. The upper horizontal axis shows the integration time re- 
quired to achieve a SNR = 5 for a star with spectral density 
r/v = (5 ■ 10" 5 photons ■ m" 2 ■ s" 1 ■ Hz" 1 ) ■ 2.5" v at | -y 1 2 = 0.5 
using a telescope of area A = 100 m 2 , a photodetector quantum 
efficiency of r) = 0.25, and a signal bandwidth of A/ = 100 MHz. 
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oped in Appendix [C] Additionally, we find that the sensi- 
tivity degradation departs from our model for large values 
of excess noise, also as expected. Finally, we see that when 
the contamination of starlight by the night sky background 
becomes important, the SNR degrades so that increasing 
the observation time no longer offers benefits to the mea- 
surements. Again we find that the SNR degrades according 
to our simple model developed in Appendix [^J These tests 
and comparisons of simulation results and simple models 
were made while isolating specific instrumental aspects one 
at a time. The good agreements obtained in all cases gives us 
confidence that the simulations may be used to characterize 
the performances of realistic instruments via their detailed 
modeling. This is important to gain an understanding of the 
sensitivity of existing and planned instruments and also this 
may be used to develop data correction algorithms so as to 
alleviate instrumental effects impacting the signals. 

Many other effects that will be encountered during real 
SII measurements can also be investigated with the sim- 
ulations. For example, the simulations can be used to in- 
vestigate the effects of inaccuracies in the time alignment 
of the signals. Furthe rmore, ACT optics ar e generally of 
Davies-Cotton design dDavies fc CottorJ[i957h . which is not 
isochronous. In previous studies, this was approximately ac- 
counted for as a signal bandwidth limitation while we could 
now simulate non- isochronous effects in detail. 

A possible implementation approach to SII consists in 
digitizing the individual telescope signals. The required fast 
digitization rate, exceeding 100 mega-samples per second, 
implies that huge volumes of data must be handled. How- 
ever, this provides a lot of flexibility as the correlation is 
obtained off-line at the time of data analysis. The simula- 
tion could also incorporate a model of the digitizer. 

The simulation algorithm presented in this paper can 
in fact be combined with any details of the envisioned in- 
strumentation to develop and optimize data analysis and 
characterize the performance of any stellar intensity inter- 
ferometer. 
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APPENDIX A: PARTIALLY CORRELATED 
POISSON STATISTICS 

The total number of photons (i.e. from the star and from 
the background light) incident on a telescope i in a time St 
is given by 

IlT, = P(n* ; (l - k) + piK + fiBGi) 

= P(fi n (l- K) + (hK + /3n n ) (Al) 

where P(m) represents a Poisson distribution of mean and 
variance m, n i , i is the mean number of photons emitted by 
the star and collected by telescope i, and nsd '■= fin i , i with 
/3 being a dimensionless factor which provides a measure of 
the amount of background light in terms of the amount of 
light from the star. 

The mean number of photons incident on a telescope 
is the sum of the photons from the star and from the night 
sky: 

n Tt = n n (1 + P) (A2) 

where n« = $* A St. 

From equations IA1I and IA2I the signal from telescope i 

is 

st { = riT i — UTf = P(ft*i(l — k) + fiiK + 0n*i) — n i , i (1 + /3) 
so that (sTj) = 0. 
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The correlation g 2 between the two telescopes signals 



g 2 _ {sn-ST?) 



Considering just the numerator: 



(A3) 



{s Tl ■ st 2 ) = ((n Tl - n Tl ) ■ (n T2 - nr 2 )) 

= ((P(n*! (1 - k) + /Zi« + /3n* J - fi*! (1 + /?)) 

•(P(n* a (1 - k) + /2 2 k + /3n* 2 ) - n* 2 (1 + /?))) 

Using Raikov's Theorem, which states that for a ran- 
dom variable Xi, if 

i, = P(rm) then ^ y Xj = P (^T m^j , (A4) 

then (sT t ■ st 2 ) can be rewritten: 

{s Ti -st 2 ) = ((P(n* 1 (l-«)) + P(/2i«)+P(/?n* 1 ) 
-^(1 + /?)) 

•(P(fi* 2 (1 - «)) + Pfa K ) + P(f3n* 2 ) 
n* 2 (l +/?))> 

Expanding and proceeding with the calculation by consid- 
ering each term separately, it is found that: 



(s Tl ■ ST 2 ) 



K (fii ■ 1x2) — n i , 1 n ir2 k 



L ) ■ 02-™* 2 )> 



(A5) 



recalling the fact that (pi) = n i , i . 

The simulation correlation is found by substituting 
equation IA5I into equation IA3I 



(an ■ st 2 ) 



^ 2 ((Mi 



(M2 



(A6) 



n 

then 



(1 + /3) 2 (1 + /3) 2 

Using the definition of the squared degree of coherence, 

2 _ - n*i) ■ (ga - n» 2 )) 



2 

5 = 



k 2 |7| 2 



(1 + /?) 2 

and solving for the quantity of interest, I7I 2 : 
, _ (1 + P? 2 



171 = 



(A7) 



(A8) 



(A9) 



where k := vAff- 

Since only the starlight may present correlation, the 
signal-to-noise expression used in Figure [5] which accounts 
for stray light contamination is derived from Equation [2] 
with 77 = 77* where 77* is the spectral density of the star, and 
Equation [3] with 77 = 77* (1 + 0) so that: 



s 

N= Aa 



(A10) 



APPENDIX B: EFFECT OF DETECTOR AREA 
ON THE COMPLEX DEGREE OF COHERENCE 

The effect of large telescope area may be quite significant 
in practice, so it is of interest to be able to predict how it 
affects and degrades the signal. To do this, first we set the 



origin of the coordinate system at the light source. Points 
on the source are labeled by positions x' with respect to the 
origin. The vector radii to the same points from a far away 
observer at position x are denoted: 



(Bl) 



For a point-like detector, the observed amplitude at x is 

A{x , t) ^(fo — uf+^(f,i))j2 / 



E(x,t) 



(B2) 



where <j>{r, t) is a random phase caused by atmospheric tur- 
bulence among other factors. The following approximation 
can be made: 

kl « N - * • |J , (B3) 
so that 

E(x,t) = ^r f A(^i) e H fc3 ''A-^+^)}dV. (B4) 

M J 

If the detector has a finite area, then the amplitude at posi- 
tion a; is a superposition of amplitudes at positions x — Xd, 
where Xd are points in the detector with respect to posi- 
tion x. The random phase can be expressed as a function 
of detector coordinates as <f>(xd, i). Now the superposition of 
amplitudes is expressed as a convolution with the detector 
area, i.e. 



E{x,t) 



C f A(^t)eW &Z ''^- ut+ ^'<>}dVd 2 



d x d Xd- 



To calculate the time averaged correlation between detectors 
i and j, denoted as (E(xi)E*(xj)), note that 

A(x, t)A(x", tyw***.*)-^*.*))^ 

= I(x')5(x — x")5(x d i — Xdj), 

where I(x) is the light intensity at point x . This is because 
separate points on the source are not correlated over large 
distances. The phase is also not correlated between separate 
points, that is, 4>(xdi,t) — <f>(xdj,t) will only be zero when 
Xdi = Xdj\ otherwise it will have a time variation which 
results in L i W s di,t)-<l>{S ii ,t))\ = Q when x ^ ^ Xdj _ 

Now defining Zi = Xi — Xdi, the time averaged correlation is 



(E{x l )E*{x ] )) = C J(f)e 



d x d Xdj- 

where C is a constant. When \x'\ <C \xj\, then the angle 9 
can be defined as 

|Sr, (B5) 
The correlation can now be expressed as: 

{E(£i)E*(i s )) = [ I(8)e-' k5 - (z >-^ d 2 9d 2 x d] (B6) 

I(zi — Sj) d 2 x d j, (B7) 



where I(zi — Zj) is the Fourier transform of the radiance 
distribution of the star, which goes from angular space to 
detector separation space. Now the quantity measured in 
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intensity interferometry is 



V\E(x^\E(x,)\ 
1 



(B8) 



I{Xi)I(Xj 



\I(zi — Zj)\ 2 d 2 x d id 2 x d j- 



Therefore, the effect of having finite sized telescopes is to 
replace the magnitude of mutual degree of coherence \y\ 2 by 
its successive convolutions with each of the telescope light 
collection area shapes. 



APPENDIX C: EXCESS NOISE 

Ignoring the effects of stray light, the number of photons 
during one time-step St in channel i is 



P(n 4 (l — «) + HiK). 



(CI) 



where P(m) represents a Poisson distribution of mean and 
variance m. 

The correlation is a result of the non-Poisson term, fn. 
The number of photons incident on channel i can also be 
written as the sum of a Poisson term and a Binomial term 
as follows: 



m = P(m(l - x)) + B(x)nj 



(C2) 



where \ ls the probability that when a photon arrives in 
channel i there is a correlated photon that also arrives in 
channel j, (rii) — fn and for simplicity we set rii £ {0, 1}, 
restricting ourselves to a situation in which (rii) <C 1. 

The excess noise is introduced as a Gaussian variation 
in the single photon response pulse. Our model simply mul- 
tiplies the amplitude by a random Gauss variable of mean 1 
and standard deviation a: 



= (P(n(l - x)) G(l, a) + B( X ) n 2 G(I, a) - n) 
•(P(n(l - X )) G(l, a) + B( X ) m G(l, a) - n). 

The effect of the excess noise on a correlated signal should be 
equivalent to the effect on an uncorrelated signal, i.e. x — 0- 
Applying that, the correlation can be rewritten as: 

si-s 2 = (P{n) G(l, a) - n) ■ (P(n) G(l, <r) - n) (C8) 
= P(n)G(l,a)P(n)G{l,a)+n 2 
-P(fi) G(l, a) n - P(n) G(l, a) n. 

Note: Suppose X and Y ar e random independent vari- 
ables then (lAhn fc Fesslerll2003l ) 

Var(XY) = E 2 (X) Var(Y) + E 2 (Y) Var(X) 
+Var(X) Var(Y) 

The variance of each term can be calculated individually, 
resulting in 

Var( Sl -s 2 ) = o i (n i + 2n 3 + n 2 ) + a 2 (2n 4 + An +2n 2 ) 
+ (2n 3 + n 2 ) + 2(a 2 (n i + n 3 ) + n 3 ) 
= cr 4 (fi 4 + 2n 3 +fi 2 ) + (j 2 (4n 4 + 6n 3 + 2fi 2 ) 



Then, 

Var(\j\ 2 ) 

And 

StdevQyf) 



+(4n 3 +n 2 ) 



Var(si ■ S2) 



-^[a (n +2n + n ) 

+a 2 (4n 4 + 6n + 2n 2 ) + (4n 3 + n 



(C9) 



(CIO) 



(C3) 



The standard deviation in equation lClOl is the standard 
deviation of just one time-step, while the term of interest is 
the standard deviation of the entire measurement, or the 
standard deviation of the mean: 



where the signal Si is made to be AC coupled by subtracting 
the mean number of photons n from the sum of individual 
photon signals (each of mean 1.0) so that the mean of the 
signal is 0. 

Since the number of photons incident on a telescope per 
time-step is considered to be small (n; G {0, 1}), the signal 
can be rewritten as 



rii G(l, a) — n. 



(C4) 



The quantity of interest is the effect that the Gaussian 
factor has on the sensitivity of the measurement, so standard 
deviation of the following term will be calculated: 



171 (m> <n 2 ) 
with 

s^(| 7 | 2 ) = stdev({si ■ S2)) 



(C5) 



(C6) 



ni ri2 

To lighten notation we assume that the two telescopes are 
identical and receive the same amount of light from the star 
(i.e. n\ =fi2 =n). 

Developing the product of the signals, 



Stdev{\i\- 



1 



Nn 2 



[a 4 (n 4 + 2ft 3 + n 2 ) 



(Cll) 



+cr 2 (4n 4 + 6n 3 + 2n 2 ) + (4n 3 + n 2 ) 



where N = is the number of time-steps taken in the 
simulation. 

Note that in the calculation described above, the Gaus- 
sian variable is truncated at zero to avoid multiplying the 
pulse by a negative number. This becomes more significant 
for large values of excess noise a while it is negligible when 
o < 1. 

This paper has been typeset from a TgX/ DTjtX file prepared 
by the author. 



(m G(l, a)-fi)- (n 2 G(l, o) - nj 



(C7) 



